"""
Week 9: VC维度分析实现
VC Dimension Analysis Implementation
"""

import numpy as np
import matplotlib.pyplot as plt
from itertools import combinations, product
import matplotlib
matplotlib.rcParams['font.sans-serif'] = ['SimHei']
matplotlib.rcParams['axes.unicode_minus'] = False

class VCDimensionAnalyzer:
    """VC维度分析器"""
    
    def __init__(self):
        self.name = "VC Dimension Analyzer"
    
    def generate_all_labelings(self, n: int) -> list[list[int]]:
        """生成n个样本的所有可能标签组合"""
        return [list(labeling) for labeling in product([0, 1], repeat=n)]
    
    def linear_classifier_2d(self, points: np.ndarray, w: np.ndarray, b: float) -> np.ndarray:
        """2D线性分类器"""
        return (points @ w + b > 0).astype(int)
    
    def can_shatter_linear_2d(self, points: np.ndarray, verbose: bool = False) -> tuple[bool, int]:
        """
        检查线性分类器是否能打散给定点集
        """
        n_points = len(points)
        all_labelings = self.generate_all_labelings(n_points)
        
        realizable_count = 0
        
        for target_labels in all_labelings:
            # 尝试找到能实现这个标签的线性分类器
            if self.find_linear_separator_2d(points, target_labels):
                realizable_count += 1
                if verbose:
                    print(f"可实现标签: {target_labels}")
        
        total_labelings = 2 ** n_points
        can_shatter = (realizable_count == total_labelings)
        
        if verbose:
            print(f"总标签数: {total_labelings}")
            print(f"可实现标签数: {realizable_count}")
            print(f"能否打散: {can_shatter}")
        
        return can_shatter, realizable_count
    
    def find_linear_separator_2d(self, points: np.ndarray, labels: list[int]) -> bool:
        """
        使用线性规划方法寻找线性分离器
        简化实现：使用感知机算法
        """
        X = np.column_stack([points, np.ones(len(points))])  # 添加偏置项
        y = np.array(labels)
        y[y == 0] = -1  # 转换为 {-1, +1}
        
        # 感知机算法
        w = np.zeros(X.shape[1])
        max_iterations = 1000
        
        for _ in range(max_iterations):
            misclassified = False
            for i in range(len(X)):
                if y[i] * (X[i] @ w) <= 0:
                    w += y[i] * X[i]
                    misclassified = True
            
            if not misclassified:
                return True
        
        return False
    
    def demonstrate_vc_dimension_linear(self):
        """演示线性分类器的VC维度"""
        print("=== 线性分类器VC维度演示 ===\n")
        
        # 1. 测试不同数量的点
        max_points = 5
        results = {}
        
        for n_points in range(1, max_points + 1):
            print(f"测试 {n_points} 个点:")
            
            # 生成随机点
            np.random.seed(42 + n_points)
            points = np.random.randn(n_points, 2)
            
            can_shatter, realizable_count = self.can_shatter_linear_2d(points, verbose=True)
            
            results[n_points] = {
                'points': points,
                'can_shatter': can_shatter,
                'realizable_count': realizable_count,
                'total_labelings': 2 ** n_points
            }
            print()
        
        # 可视化结果
        self.visualize_vc_results(results)
        
        return results
    
    def visualize_vc_results(self, results: dict):
        """可视化VC维度结果"""
        fig, axes = plt.subplots(2, 3, figsize=(18, 12))
        
        # 1-5个点的可视化
        for i, (n_points, result) in enumerate(results.items()):
            if i >= 5:
                break
                
            row = i // 3
            col = i % 3
            ax = axes[row, col]
            
            points = result['points']
            can_shatter = result['can_shatter']
            
            # 绘制点
            ax.scatter(points[:, 0], points[:, 1], s=100, c='red', 
                      edgecolors='black', linewidth=2)
            
            # 标注点的索引
            for j, point in enumerate(points):
                ax.annotate(f'{j}', (point[0], point[1]), 
                           xytext=(5, 5), textcoords='offset points')
            
            ax.set_title(f'{n_points}个点 - {"能打散" if can_shatter else "不能打散"}')
            ax.grid(True, alpha=0.3)
            ax.set_aspect('equal')
        
        # 第6个子图：增长函数
        ax = axes[1, 2]
        
        n_points_list = list(results.keys())
        realizable_counts = [results[n]['realizable_count'] for n in n_points_list]
        total_labelings = [results[n]['total_labelings'] for n in n_points_list]
        
        ax.plot(n_points_list, realizable_counts, 'bo-', linewidth=2, 
               markersize=8, label='实际可实现数')
        ax.plot(n_points_list, total_labelings, 'ro--', linewidth=2, 
               markersize=8, label='总标签数 (2^n)')
        
        ax.set_xlabel('点数 n')
        ax.set_ylabel('标签数')
        ax.set_title('增长函数')
        ax.legend()
        ax.grid(True, alpha=0.3)
        ax.set_yscale('log')
        
        plt.tight_layout()
        plt.show()

class GrowthFunctionAnalyzer:
    """增长函数分析器"""
    
    def __init__(self):
        self.name = "Growth Function Analyzer"
    
    def theoretical_growth_function(self, n: int, vc_dim: int) -> int:
        """理论增长函数上界"""
        if n <= vc_dim:
            return 2 ** n
        else:
            # Sauer引理: m_H(n) <= sum_{i=0}^{vc_dim} C(n,i)
            total = 0
            for i in range(vc_dim + 1):
                total += self.binomial_coefficient(n, i)
            return total
    
    def binomial_coefficient(self, n: int, k: int) -> int:
        """计算二项式系数 C(n,k)"""
        if k > n or k < 0:
            return 0
        if k == 0 or k == n:
            return 1
        
        # 使用递推公式避免大数计算
        result = 1
        for i in range(min(k, n - k)):
            result = result * (n - i) // (i + 1)
        
        return result
    
    def compare_growth_functions(self, max_n: int = 20):
        """比较不同VC维度的增长函数"""
        n_values = range(1, max_n + 1)
        vc_dims = [1, 2, 3, 5, 10]
        
        plt.figure(figsize=(15, 10))
        
        # 增长函数比较
        plt.subplot(2, 2, 1)
        
        for vc_dim in vc_dims:
            growth_values = [self.theoretical_growth_function(n, vc_dim) for n in n_values]
            plt.plot(n_values, growth_values, 'o-', label=f'VC维={vc_dim}', linewidth=2)
        
        # 添加指数函数作为参考
        exponential = [2**n for n in n_values]
        plt.plot(n_values, exponential, 'k--', alpha=0.7, label='2^n')
        
        plt.xlabel('样本数 n')
        plt.ylabel('增长函数值')
        plt.title('不同VC维度的增长函数')
        plt.legend()
        plt.grid(True, alpha=0.3)
        plt.yscale('log')
        
        # 增长函数比率
        plt.subplot(2, 2, 2)
        
        for vc_dim in vc_dims[:-1]:  # 排除最大的VC维度
            growth_values = [self.theoretical_growth_function(n, vc_dim) for n in n_values]
            exponential = [2**n for n in n_values]
            ratios = [g/e for g, e in zip(growth_values, exponential)]
            plt.plot(n_values, ratios, 'o-', label=f'VC维={vc_dim}', linewidth=2)
        
        plt.xlabel('样本数 n')
        plt.ylabel('增长函数 / 2^n')
        plt.title('增长函数与指数函数的比率')
        plt.legend()
        plt.grid(True, alpha=0.3)
        
        # VC维度与样本复杂度
        plt.subplot(2, 2, 3)
        
        # 根据VC理论，样本复杂度大约为 O(vc_dim * log(1/delta) / epsilon^2)
        epsilon = 0.1
        delta = 0.1
        
        sample_complexities = []
        for vc_dim in range(1, 21):
            # 简化的样本复杂度公式
            sample_complexity = (vc_dim * np.log(1/delta)) / (epsilon**2)
            sample_complexities.append(sample_complexity)
        
        plt.plot(range(1, 21), sample_complexities, 'bo-', linewidth=2)
        plt.xlabel('VC维度')
        plt.ylabel('样本复杂度')
        plt.title(f'样本复杂度 vs VC维度 (ε={epsilon}, δ={delta})')
        plt.grid(True, alpha=0.3)
        
        # 泛化界
        plt.subplot(2, 2, 4)
        
        # 泛化界: sqrt(vc_dim * log(n) / n)
        n_samples = np.logspace(1, 4, 50)
        
        for vc_dim in [1, 5, 10, 20]:
            generalization_bounds = np.sqrt(vc_dim * np.log(n_samples) / n_samples)
            plt.loglog(n_samples, generalization_bounds, label=f'VC维={vc_dim}', linewidth=2)
        
        plt.xlabel('样本数 n')
        plt.ylabel('泛化界')
        plt.title('VC泛化界')
        plt.legend()
        plt.grid(True, alpha=0.3)
        
        plt.tight_layout()
        plt.show()

class ShatteringExperiment:
    """打散实验类"""
    
    def __init__(self):
        self.name = "Shattering Experiment"
    
    def interval_classifier_1d(self, points: np.ndarray, interval: tuple[float, float]) -> np.ndarray:
        """1D区间分类器"""
        a, b = interval
        return ((points >= a) & (points <= b)).astype(int)
    
    def can_shatter_intervals_1d(self, points: np.ndarray) -> tuple[bool, int]:
        """检查区间分类器是否能打散1D点集"""
        n_points = len(points)
        points_sorted = np.sort(points)
        all_labelings = self.generate_all_labelings(n_points)
        
        realizable_count = 0
        
        for target_labels in all_labelings:
            if self.find_interval_separator_1d(points, target_labels):
                realizable_count += 1
        
        total_labelings = 2 ** n_points
        can_shatter = (realizable_count == total_labelings)
        
        return can_shatter, realizable_count
    
    def generate_all_labelings(self, n: int) -> list[list[int]]:
        """生成所有标签组合"""
        return [list(labeling) for labeling in product([0, 1], repeat=n)]
    
    def find_interval_separator_1d(self, points: np.ndarray, labels: list[int]) -> bool:
        """寻找区间分离器"""
        # 对于区间分类器，正例必须形成连续区间
        positive_indices = [i for i, label in enumerate(labels) if label == 1]
        
        if len(positive_indices) == 0:
            return True  # 空区间
        
        # 检查正例是否连续
        points_with_indices = [(points[i], i) for i in range(len(points))]
        points_with_indices.sort()  # 按点的值排序
        
        sorted_labels = [labels[idx] for _, idx in points_with_indices]
        
        # 检查正例是否形成连续段
        in_positive_region = False
        positive_started = False
        
        for label in sorted_labels:
            if label == 1:
                if not in_positive_region and positive_started:
                    return False  # 正例不连续
                in_positive_region = True
                positive_started = True
            else:
                if in_positive_region:
                    in_positive_region = False
        
        return True
    
    def demonstrate_interval_shattering(self):
        """演示区间分类器的打散能力"""
        print("=== 区间分类器打散实验 ===\n")
        
        results = {}
        
        for n_points in range(1, 6):
            print(f"测试 {n_points} 个点:")
            
            # 生成排序的点
            points = np.linspace(0, 1, n_points)
            
            can_shatter, realizable_count = self.can_shatter_intervals_1d(points)
            
            results[n_points] = {
                'points': points,
                'can_shatter': can_shatter,
                'realizable_count': realizable_count,
                'total_labelings': 2 ** n_points
            }
            
            print(f"  点集: {points}")
            print(f"  能否打散: {can_shatter}")
            print(f"  可实现标签数: {realizable_count}/{2**n_points}")
            print()
        
        # 可视化
        self.visualize_interval_shattering(results)
        
        return results
    
    def visualize_interval_shattering(self, results: dict):
        """可视化区间打散结果"""
        fig, axes = plt.subplots(2, 3, figsize=(18, 8))
        
        # 前5个子图显示不同点数的情况
        for i, (n_points, result) in enumerate(results.items()):
            if i >= 5:
                break
                
            row = i // 3
            col = i % 3
            ax = axes[row, col]
            
            points = result['points']
            can_shatter = result['can_shatter']
            
            # 绘制点
            ax.scatter(points, np.zeros_like(points), s=100, c='red', 
                      edgecolors='black', linewidth=2)
            
            # 标注点
            for j, point in enumerate(points):
                ax.annotate(f'{j}', (point, 0), xytext=(0, 10), 
                           textcoords='offset points', ha='center')
            
            ax.set_xlim(-0.1, 1.1)
            ax.set_ylim(-0.5, 0.5)
            ax.set_title(f'{n_points}个点 - {"能打散" if can_shatter else "不能打散"}')
            ax.grid(True, alpha=0.3)
        
        # 第6个子图：增长函数比较
        ax = axes[1, 2]
        
        n_points_list = list(results.keys())
        realizable_counts = [results[n]['realizable_count'] for n in n_points_list]
        total_labelings = [results[n]['total_labelings'] for n in n_points_list]
        
        # 理论值：区间分类器的VC维度是2
        theoretical_values = [min(2**n, n+1) for n in n_points_list]
        
        ax.plot(n_points_list, realizable_counts, 'bo-', linewidth=2, 
               markersize=8, label='实际可实现数')
        ax.plot(n_points_list, theoretical_values, 'go-', linewidth=2, 
               markersize=8, label='理论值 (VC维=2)')
        ax.plot(n_points_list, total_labelings, 'ro--', linewidth=2, 
               markersize=8, label='总标签数 (2^n)')
        
        ax.set_xlabel('点数 n')
        ax.set_ylabel('标签数')
        ax.set_title('区间分类器增长函数')
        ax.legend()
        ax.grid(True, alpha=0.3)
        
        plt.tight_layout()
        plt.show()

def demonstrate_vc_dimension():
    """演示VC维度概念"""
    print("=== VC维度概念演示 ===\n")
    
    # 1. 线性分类器VC维度
    vc_analyzer = VCDimensionAnalyzer()
    linear_results = vc_analyzer.demonstrate_vc_dimension_linear()
    
    # 2. 增长函数分析
    print("=== 增长函数分析 ===\n")
    growth_analyzer = GrowthFunctionAnalyzer()
    growth_analyzer.compare_growth_functions(max_n=15)
    
    # 3. 区间分类器打散实验
    print("=== 区间分类器实验 ===\n")
    shattering_exp = ShatteringExperiment()
    interval_results = shattering_exp.demonstrate_interval_shattering()
    
    # 总结
    print("=== VC维度理论总结 ===")
    print("1. 线性分类器在2D空间的VC维度是3")
    print("2. 区间分类器在1D空间的VC维度是2")
    print("3. VC维度决定了函数类的表达能力")
    print("4. 更高的VC维度意味着更强的拟合能力，但也可能导致过拟合")
    print("5. 样本复杂度与VC维度成正比")

if __name__ == "__main__":
    demonstrate_vc_dimension()